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Diverse complex dynamical systems are known to exhibit abrupt regime shifts at bifurcation 
points of the saddle-node type. The dynamics of most of these systems, however, have a stochastic 
component resulting in noise driven regime shifts even if the system is away from the bifurcation 
points. In this letter, we propose a new quantitative measure, namely, the propensity transition 
point as an indicator of stochastic regime shifts. The concepts and the methodology are illustrated 
for the May model, a well-known model in ecology. The general applicability and usefulness of the 
method for the analysis of regime shifts is further pointed out. 
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INTRODUCTION 



Dynamical systems are known to undergo sudden 
regime shifts at critical parameter values, termed the bi- 
furcation points, with the different regimes defined by 
distinct sets of attractors of the dynamics Recently, 
a large number of studies have been devoted to the inves- 
tigation of regime shifts in complex dynamical systems 
and processes ranging from ecosystems, climate change, 
population collapse and financial markets to epilep- 
tic seizures and asthma attacks 0, 0] ■ Examples of sud- 
den regime shifts include the collapse of vegetation un- 
der semi-arid conditions, the transition from a clear to a 
turbid lake, catastrophic shifts in fish or wildlife popu- 
lations [2t3 , transitions from one stable gene expression 
state to another in natural and synthetic genetic circuits 
SHlOj and the sudden deterioration of complex diseases 
The regime shifts, in most of the cases studied, are 
brought about by abrupt transitions from bistability to 
monostability via the saddle- node bifurcation (T], Q ■ The 
state of a dynamical system at a specific time t is defined 
in terms of the magnitudes of one or more key variables 
at time t. In the steady state, the rates of change of the 
variables are zero. In the following, we focus our atten- 
tion on a single variable dynamical system. In the case 
of monostability, there is only one single stable steady 
state. The region of bistability is distinguished by the 
existence of two stable steady states separated by an un- 
stable steady state. The two stable steady states corre- 
spond to low and high values of the variable, designated 
as the L and H states respectively. The bistable region 
separates two regions of monostability with the L and 
H states being the respective stable steady states. At a 
saddle node bifurcation point, one of the stable steady 
states merges with the unstable steady state leading to 
the disappearance of the physical solutions beyond the 
bifurcation point. Figure [1] illustrates the bifurcation in 
the case of the May model [l^l, a well-known model in 
ecology. In the original model, a population of herbivores 
at a constant density is considered, the population being 
sustained by vegetation of biomass K. The rate equation 



describing the dynamics of the biomass is given by 
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where x represents the biomass. The first term on 
the right hand side describes the logistic growth of the 
biomass at the rate r, K being the carrying capacity with 
the growth rate becoming zero when x = K . The second 
term corresponds to the rate of loss of the biomass due to 
the consumption of the biomass by the herbivores. The 
parameter c is the maximum loss rate of the biomass or 
equivalently the maximum grazing rate. The loss rate 
saturates when x is much larger than a characteristic 
value xq. The May model has been applied to various 
other ecological problems, e.g., the exploitation of fish 
populations [l^ , dynamics of spruce budworms [l^ and 
harvesting of macrophytes [l^ . Figure[T]shows the steady 
state values of the biomass x = O) versus the maxi- 
mum grazing rate parameter c. The other parameter val- 
ues are r = \, K = \Q and = 1 in appropriate units. 
The solid lines represent the stable steady states sepa- 
rated by a branch (dot-dashed line) of unstable steady 
states. The parameters chosen fall in a region of param- 
eter space exhibiting bistability. The lower (ci) and up- 
per (C2) bifurcation points separate a region of bistability 
from two regions of monostability. At the upper bifurca- 
tion point C2 defined by a critical value of c, a saddle node 
bifurcation take place and there is an abrupt regime shift 
from the H (high value of biomass) to the L (low value of 
biomass) state. The transition to the collapsed biomass 
state is not reversible with the reverse transition from the 
L to the H state occurring at the lower bifurcation point 
ci. This is the phenomenon of hysteresis, a characteristic 
feature of regime shifts of the type shown in Figure [T] 
At ci, a saddle node bifurcation occurs with the L state 
losing stability. 

In the case of deterministic time evolution of a dynam- 
ical system, there are two ways in which regime shifts 
occur: (i) at the bifurcation points and (ii) on apply- 
ing large perturbations in the region of bistability. The 
unstable steady state sets a threshold for switching tran- 
sitions between the L and H states. Consider the system 
to be initially in the H state. A sufficiently large pertur- 
bation reduces the magnitude of the dynamical variable 
below the threshold, thereby bringing about a transition 
to the L state. The time evolution in complex dynamical 
systems exhibiting reg ime shifts is, in general, stochas- 
tic in nature 0, [l3 - ll6j due to the probabilistic nature of 
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Figure 1. Biomass x versus the maximum grazing rate pa- 
rameter c in the steady state of the May model tlie dynamics 
of which are described by Eq. (1). The solid lines represent 
the stable steady states and the dot-dashed line the branch of 
unstable steady states. The points c = ci and c = C2 repre- 
sent the lower and upper bifurcation points respectively. The 
parameter values used are r = 1, K = 10 and xo = 1. The in- 
set shows the plot of the steady state probability distribution 
Pst(x) versus x in the presence of additive noise of strength 
di = .25. The maximum grazing rate parameter c = 2.3. The 
extrema points correspond to the stable steady states of the 
deterministic dynamics, with x\ and 13 representing the sta- 
ble steady states which are separated by an unstable steady 
state X2- 

the processes associated with the dynamics. Stochastic- 
ity gives rise to fluctuations in the magnitude of the dy- 
namical variable (say, the biomass) which, if sufficiently 
strong, cause transitions between the L and H states 
away from the bifurcation points. The occasional switch 
to an alternative regime in the presence of stochasticity is 
known as "flickering" Q and is responsible for the appear- 
ance of a bimodal steady state probability distribution 
in the region of "bistability" . Noise-induced excursions 
from the stable steady states result in a broadening of 
the distribution around the steady states. The steady 
state distribution is achieved when the probability of the 
H to the L state transition is the same as that of the 
reverse transition. The bifurcation theory of stochastic 
dynamical systems is not as well-formulated as that in 
the case of deterministic dynamics. Kepler and Elston 
[13] have defined stochastic bifurcation in terms of the 
number of critical (singular) points in the steady state 
probabiHty density function. Song et. al. [3 have in- 
vestigated the stochastic bifurcation structure of cellular 
networks through specifying four quantities as a function 
of the bifurcation parameter: (i) the number of distinct 
subpopulations (or peaks in the steady state probability 
distributions), (ii) the locations of the peaks in terms of 
the measurable variable x, (in) the variability in x for 
each subpopulation and (iv) the fractions of the whole 
population represented by the subpopulations. The total 
probability distribution is expressed as a mixture of com- 
ponent distributions, one for each subpopulation. The 
term subpopulation may be used in a generalized sense. 
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Figure 2. The equipartition point 2:1/2 (Eq. (14)) versus the 
parameter c for three different strengths of the additive noise, 
di = .15 (solid line), di = .20 (dashed line) and di = .25 
(dot-dashed line). The other parameter values are r = 1, 

K = IG and xq = 1. The inset shows the plots of versus 
c for the three noise strengths. The peak positions define the 
propensity transition point c*. 



the number of subpopulations being equal to the number 
of distinct peaks in the probability distribution. For ex- 
ample, in the region of bistability (Figure [T|), one has two 
subpopulations, L and H and the corresponding proba- 
bility distributions arise due to noise-induced broade ning 
of the steady state levels. Guttal and Jayaprakash [l^l 
have investigated the impact of noise on bistable ecolog- 
ical systems. They have shown that the region of bista- 
bility is diminished in the presence of small amounts of 
noise whereas for noise beyond a critical strength, bista- 
bility vanishes altogether. In the latter case, the dynam- 
ical system can undergo abrupt regime shifts frequently. 
In this Letter, we propose a new framework for the quan- 
titative characterization of noise-induced regime shifts in 
bistable systems. The methodology has a clear physi- 
cal interpretation and is easy to implement. We consider 
the May model to illustrate the methodology which is of 
general applicability. 



II. EQUIPARTITION AND PROPENSITY 
TRANSITION POINTS 

In the case of stochastic dynamics, the one-variable 
Langevin equation (LE) including both additive and mul- 
tiplicative noise terms is given by 

^ = /(x)+5(x)e(0+r(i) (2) 

where T{t) (additive noise) and £{t) (multiplicative noise) 
are Gaussian white noises with zero mean and correla- 
tions: 

{T{t)T{t')) = 2diS(t - t') 
{e{t)e{t')) = 2d25(t - t') (3) 

In Eq. ([3]), di and ^2 denote the strengths of the addi- 
tive and multiplicative noises respectively, we further as- 
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sume that there is no correlation between the two types 
of noise. The first term in the right hand side of Eq. ([5]) 
describes the deterministic dynamics, e.g., for the May 
model f{x) is given by the right hand side expression in 
Eq. ([T]). The Fokker-Planck (FP) equation, a rate equa- 



tion for the probability distribution P(x,t), is obtained 
mr the L 



jbabiiity 
from the LE as [13, 



dt 



— t)] + — [B{x)P{x, t)] (4) 



dx 



where 



and 



A{x)= f{x)+d2g{x)g'{x) 



B{x)^di+d2[g{x)f 



(5) 



(6) 



The steady state probability distribution (SSPD) is given 
by ^21.2^ 



Pst{x) 



N 



■ exp 



Mxl 

B{x) 



dx 



N 



{d2[ff(^)]2+di}^ 

f{x')dx' 



X exp 



d2[9{x')f + di 



(7) 



where N is the normalization constant. Equation Q can 
be rewritten in the form 



where 



(I)f{.x) = - In [d2[g{x)f + di] 



f{v)dy 



(8) 



d2 [9[y)Y + d^ 



(9) 



defines the "stochastic potential" of the dynamics. The 
inset of Figure [T] shows the SSPD, Pst{x), versus x for the 
May model in the presence of only additive noise {g{x) = 

0, d2 = 0) for the parameter values r = 1, K = 10, xq = 

1, c = 2.3 and di = 0.25. In this case, 4)f{x) reduces to 

^p{x) = ^lnd,-j-J f{y)dy (10) 
The deterministic potential (x) is given by 

fiy)dy (11) 



Comparing Eqs. PH)) and (|lip one finds that the cxtrema 
points of (f>F{x) and (poix) are identical with the minima 
(maximum) corresponding to stable (unstable) steady 
states. In the SSPD Pst{x), the stable steady states cor- 
respond to the maxima points xi and x^ whereas the min- 
imum point X2 represents the unstable steady state. The 




Figure 3. The mean first passage times tl (solid line) and th 
(dashed line) versus the parameter c. The other parameter 
values are r = 1, K = 10, xo = 1 and di = .15. The inter- 
section point Ct {tl = Th) of the two curves is approximately 
equal to c*, the propensity transition point (c-r = 2.360, 
c* = 2.364). 



cumulative distribution function (CDF) Pc{x) is given by 



Pc {x') 



P{x)dx 



(12) 



where < a;' < oo. Let Xi/2 be the value of x which 
marks the equipartition point of the SSPD, Pst{x), i.e., 



P{x)dx 



P{x)dx 



Xl/2 



Pc (2^1/2) = 1 - Pc (2:^1/2) 



(13) 



(14) 



As the bifurcation parameter c is changed, the SSPD 
Pst{x) also changes from which the corresponding Xi/2 
can be determined using Eq. ([T^ . Figure [2] shows the 
variation of Xi/2 with c for the parameter values r = 1, 
= 10, Xo — 1 and different values of the noise strength 
di = .15 (solid line), di = .2 (dashed line) and di = .25 
(dot-dashed line). The deterministic bifurcation points 
have values ci ~ 1.788 and C2 = 2.604. As c is increased 
from a low value, one notices a transition of Xi/2 from a 
high to a low value. The transition is sharp for low values 
of the additive noise and becomes smeared as the noise 
strength di is increased. The point x = Xi/2 provides a 
quantitative measure of which subpopulation is the pre- 
ferred (dominant) subpopulation in the steady state. The 
L subpopulation corresponds to x values in the range 
< X < X2, the rest of the SSPD is associated with the 
H subpopulation. For low values of c, Xi/2 is > X2 and 
the H subpopulation is more dominant. As c increases, 
Xi/2 shifts towards X2 and the L subpopulation becomes 
the more preferred one when Xi/2 is < a;2. One can de- 
fine a propensity transition point (PTP) c* at which the 
preferred subpopulation changes its character from H to 
L. The inset of Figure [5] shows the variation of 
versus c for di = .15, di = .2 and di = .25, with the 
respective peak positions being c* = 2.364, c* = 2.371 
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(a) 



di=0.15, d2=0.0 ■ 
di=0.15, d,=0.1 




Figure 4. (a) The equipartition point x-^j^ (Eq. (O) versus c 
for different strengths of the multiplicative noise: d2 = (solid 
line), di = .1 (dashed line) with d\ = .15 in each case. The 
other parameter values are r = 1, if = 10 and a-Q = 1. (b) 
The cumulative residual entropy (CRE) (Eq. (|15|) ) versus c 
for the two different values of the multiplicative noise strength 
as in (a). The FTP c* is included for comparison. 



Pn^t) are given by 

dt 
dPnit) 
dt 



^-kLPL{t) + kHPH{t) 
= -kHPH{t) + kLPL{t) 



(16) 



where and ku are the stochastic transition rates for 
the L ^ H and H ^ L transitions respectively. In the 

dPL(t) _ ^ dP„{t) 
dt 



steady state, 



0, T- 

Phs ^ kL 
Pls kn 



and one gets 



(17) 



where '5" in the suffix indicates the steady state proba- 
bilities. In the case of weak noise, one can write [25| 



1 

TL ~ — 
kL 
1 

kn 



(18) 



where tl and th denotes the mean first passage times 
(MFPTs) for exits from the domains of the L and H 
subpopulations respectively. The exits arc brought about 
by noise and a higher value of the MFPT indicates a 
greater stability of the steady state from which the exit 
occurs. At the FTP c*, Phs = Pls which leads to the 
equality tl ^ th- For weak noise, fcff and fc^ are given 
by the approximate expressions , 

kH = ^JW^.{x2)\^'^ (a;3)e^H-3)-0H-2) 

ZTT V 



kL = i^JW^{x2)\<l>'^ (^^)e^Hxi)-0H-2) (19) 

ZTT V 

where (j)'p (X2) denotes the double derivative of (pp with 
respect to x. Figure [3] shows the plots of the MFPTs 
Tl and th (determined using Eqs. PH]) . and ([TO]) ) 
versus the bifurcation parameter c. The other parameter 
values are r = 1, A' = 10, xq = 1 and di = .15. The inter- 
section point of the two curves, Cr, at which tl = th, is 
approximately equal to c* (c^ = 2.360, c* = 2.364). The 
numerical result supports the analytic argument, using 
Eqs. (HSl) and that tl = th at the FTP c = c* . 



and c* = 2.378. For larger magnitudes of noise below 
some critical strength, the variation is more smeared but 
one can still identify a peak position c* at which 
attains its maximum value. In the presence of large ad- 
ditive noise, flickering becomes prominent and bistability 
is destroyed. In the presence of only additive noise, c* is 
equal to X2, the minimum of the SSPD, which is again 
equal to the solution for the deterministic unstable steady 
state. The cumulative residual entropy (CRE) defined as 

M 

Sc = - [ Pc{x)\nP,{x)dx (15) 

J X 

with Pc{x) as given in Eq. (jl2p . attains its maximum 
value at c = Cr which is close to c*. For di = .15, 
Cr = 2.343. Let Pl and Ph be the probabilities of 
belonging to subpopulations L and H respectively. Pl 
and PfL can be written as PL{t) ~ J^^ P{x,t)dx and 
Pnit) = /^°° P{x,t)dx. The time evolutions of PL{t) and 



We next consider the general case when additive and 
multiplicative noise terms are present in the LE, Eq. (0). 
We assume that the multiplicative noise is associated with 
the rate constant r in Eq. ([T|) with g{x) (Eq. ^) given 

by 

gix)^x[l-^) (20) 

In this case, the extrema points xi, X2 and 2:3 of the 
SSPD no longer coincide with the deterministic steady 
state values. As before, the equipartition point 2:1/2 can 
be determined as a function of c. Figure S] (a) shows 
the variation of 2:1/2 versus c for the parameter values 
r ^ I, K ^ 10 and xq = 1. Figure H(b) shows the CRE 
(Eq. (jl5p ) versus the bifurcation parameter c. The peak 
position Cr is close to the FTP c* in the presence of only 
additive noise. When both additive and multiplicative 
noise terms appear in the LE (Eq. ([2])), the difference 
between cr and c* increases. 
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III. CONCLUDING REMARKS 

The quantitative measure of regime shifts, namely, the 
PTP, proposed in this Letter is of broad applicability and 
can be computed both when the model dynamics are 
known or when only the time series data are available. 
The physical interpretation of the PTP in terms of a tip- 
ping of the balance towards one regime or the other helps 
in identifying the parameter regime to be avoided in or- 
der to forestall disastrous regime shifts. The variation of 
the cquipartition point Xi/2 versus a bifurcation param- 
eter provides advanced knowledge of the approach to the 
PTP. A number of signatures of impending regime shifts 
have been proposed so far in the scenario of the saddle- 
node bifurcation [l|, [3, [3, [3 ■ These include the criti- 



cal slowing down, rising variance and skewness of a sub- 
population probability distribution and the rising lag-1 
auto-correlation function, as a function of the bifurcation 
parameter. Some of these quantities would also provide 
characteristic signatures as the PTP is approached since 
the PTP precedes the deterministic bifurcation point. A 
comprehensive analysis of some of these issues is post- 
poned till a future publication. 
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